{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Central Limit Theorem"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Introduction\n",
    "\n",
    "This is a matplotlip example that visualizes the [**Central Limit Theorem**](https://en.wikipedia.org/wiki/Central_limit_theorem). This states that the sum of any independent random variables converges towards a normal distribution (under proper normalization). It is **not** required that the individual variables are normally distributed.\n",
    "\n",
    "Let's visualize this with a simple interactive *matplotlib* figure.\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We'll import *NumPy* and *matplotlib* and use two of the \"magic Jupyter commands\" to inline plots in a ```.svg``` vector format."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:13.939460Z",
     "start_time": "2018-10-03T13:21:13.301433Z"
    }
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "%matplotlib inline\n",
    "\n",
    "%config InlineBackend.figure_format = 'svg'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can easily generate a huge array of (uniform) random variables with *NumPy*:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:13.946898Z",
     "start_time": "2018-10-03T13:21:13.940903Z"
    }
   },
   "outputs": [],
   "source": [
    "samples = np.random.uniform(size=10)\n",
    "print(samples)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's plot such a distribution with a histogram:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:14.263141Z",
     "start_time": "2018-10-03T13:21:13.948323Z"
    },
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "plt.hist(np.random.uniform(size=100000))\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Useful parameters for histograms are:\n",
    "* *```normed```* to get a normalized output\n",
    "* *```bins```* to specify the number of bins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:14.609791Z",
     "start_time": "2018-10-03T13:21:14.265439Z"
    }
   },
   "outputs": [],
   "source": [
    "plt.hist(np.random.uniform(size=100000), bins=50, density=True)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To check the *central limit theorem* however, we need to consider the distribution of the **sums** of many random variables. This is also easy in NumPy by generating a 2d-array (or a matrix) of random values and then summing over one of the dimensions.\n",
    "Let's do the sum of $n$ random variables. For large $n$, it should converge to a normal distribution!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:14.975539Z",
     "start_time": "2018-10-03T13:21:14.612654Z"
    }
   },
   "outputs": [],
   "source": [
    "n = 10\n",
    "samples = np.random.uniform(size=(100000, n))\n",
    "sums = np.sum(samples, axis=1)\n",
    "print(sums)\n",
    "\n",
    "plt.hist(sums, bins=50, density=True)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's add a title and make the figure interactive so we can play with $n$.\n",
    "\n",
    "We import another module that enables interactivity, and then wrap the cell above into a function that has a *```@interact```* decorator."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:15.445266Z",
     "start_time": "2018-10-03T13:21:14.978753Z"
    }
   },
   "outputs": [],
   "source": [
    "from ipywidgets import interact\n",
    "\n",
    "@interact(n=(1, 20))\n",
    "def central_limit_theorem(n=1):  # Default parameter is also default slider position\n",
    "    # Look at the sum of n realizations of sampling many random numbers from any distribution (here uniform)\n",
    "    samples = np.random.uniform(size=(100000, n))\n",
    "    sums = np.sum(samples, axis=1)\n",
    "    \n",
    "    # Plot normalized histogram\n",
    "    plt.hist(sums, bins=50, density=True)\n",
    "    \n",
    "    # Set a title\n",
    "    plt.title(\"Central limit theorem\")\n",
    "    \n",
    "    # Display the figure\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To verify that the samples are indeed normally distributed, we can also overlay a true normal distribution.\n",
    "This is given by:\n",
    "\n",
    "$$f(x \\; | \\; \\mu, \\sigma^2) = \\frac{1}{\\sqrt{2\\pi\\sigma^2}} e^{-\\frac{(x - \\mu)^2}{2\\sigma^2}}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:15.458411Z",
     "start_time": "2018-10-03T13:21:15.447166Z"
    }
   },
   "outputs": [],
   "source": [
    "def normal_distr(x, mu, sigma):\n",
    "    return 1/(np.sqrt(2 * np.pi * sigma**2)) * np.exp( -(x - mu)**2 / (2 * sigma**2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:16.193994Z",
     "start_time": "2018-10-03T13:21:15.461761Z"
    },
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "from ipywidgets import interact\n",
    "\n",
    "@interact(n=(1, 20))\n",
    "def central_limit_theorem(n=1):\n",
    "    # Look at the sum of n realizations of sampling many random numbers from any distribution (here uniform)\n",
    "    samples = np.random.uniform(size=(100000, n))\n",
    "    sums = np.sum(samples, axis=1)\n",
    "\n",
    "    # Plot normalized histogram & extract the x position of the bins\n",
    "    _, x, _ = plt.hist(sums, bins=50, density=True)\n",
    "    \n",
    "    # Estimate parameters for normal distribution\n",
    "    mu = np.mean(sums)\n",
    "    sigma = np.std(sums)\n",
    "    \n",
    "    # Plot corresponding normal distribution\n",
    "    y = normal_distr(x, mu, sigma)\n",
    "    plt.plot(x, y, linewidth=2, color='red', linestyle='--')\n",
    "    \n",
    "    # Set a title\n",
    "    plt.title(\"Central limit theorem\")\n",
    "    \n",
    "    # Display the figure\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, let's add a legend by adding labels to both the *```hist```* and *```plot```* commands:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2018-10-03T13:21:16.819543Z",
     "start_time": "2018-10-03T13:21:16.198121Z"
    }
   },
   "outputs": [],
   "source": [
    "from ipywidgets import interact\n",
    "\n",
    "@interact(n=(1, 20))\n",
    "def central_limit_theorem(n=1):\n",
    "    # Look at the sum of n realizations of sampling many random numbers from any distribution (here uniform)\n",
    "    samples = np.random.uniform(size=(100000, n))\n",
    "    sums = np.sum(samples, axis=1)\n",
    "    \n",
    "    # Plot normalized histogram & extract the x position of the bins\n",
    "    _, x, _ = plt.hist(sums, bins=50, density=True, label='Approximation')\n",
    "    \n",
    "    # Estimate parameters for normal distribution\n",
    "    mu = np.mean(sums)\n",
    "    sigma = np.std(sums)\n",
    "    \n",
    "    # Plot corresponding normal distribution\n",
    "    y = normal_distr(x, mu, sigma)\n",
    "    plt.plot(x, y, linewidth=2, color='red', linestyle='--', label='Normal distribution')\n",
    "    \n",
    "    # Set a title\n",
    "    plt.title(\"Central limit theorem\")\n",
    "    \n",
    "    # Add a legend\n",
    "    plt.legend(loc='upper left')\n",
    "    \n",
    "    # Display the figure\n",
    "    plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  },
  "toc": {
   "nav_menu": {
    "height": "12px",
    "width": "252px"
   },
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": true,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": "block",
   "toc_window_display": false
  },
  "varInspector": {
   "cols": {
    "lenName": "30",
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}